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abstract 

High-resolution two-dimensional simulations were performed for the first five minutes of the evolution of a core 
collapse supernova explosion in a 15 M.q blue supergiant progenitor. The computations start shortly after bounce 
and include neutrino-matter interactions by using a light-bulb approximation for the neutrinos, and a treatment of 
the nucleosynthesis due to explosive silicon and oxygen burning. We find that newly formed iron-group elements 
are distributed throughout the inner half of the helium core by Rayleigh-Taylor instabilities at the Ni-nSi/O and 
C-nO/He interfaces, seeded by convective overturn during the early stages of the explosion. Fast moving nickel 
mushrooms with velocities up to ^ 4000 km/s are observed. This offers a natural explanation for the mixing 
required in light curve and spectral synthesis studies of Type lb explosions. A continuation of the calculations to 
later times, however, indicates that the iron velocities observed in SN 1987 A cannot be reproduced because of a 
strong deceleration of the clumps in the dense shell left behind by the shock at the He/H interface. 

Subject headings: hydrodynamics - instabilities - nucleosynthesis - shock waves - stars: supernovae 



1. INTRODUCTION 

Overwhelming observational evidence (cf. references in 



Miiller 1998) suggests that large-scale mixing processes took 
place in SN 1987 A and transported newly synthesized ^^Ni 
from its creation site near the collapsed core all the way out to 
the hydrogen envelo pe. Spectroscopic studies of SN 1987 F, SN 
1988 A, SN 199 3 J (|Spyromilio, 1994 and references therein) 
and SN 1995 V (Fassia et al., 1998) indicate that such mixing 



is probably generic in core collapse supernovae. Indeed, artifi- 
cial mixing of the radioactive ejecta within the (helium) enve- 
lope is indispensable in order to reproduce the light curves and 
spectra of Type lb e xplosions using one-dimensional (ID) hy- 
drodynamic models (jShigeyama, Nomoto, & Tsujimoto 1990 ; 



Woosley & Eastman 1997, and references therein) 

These findings have instigated theoretical work on multidi- 
mensional supernova models focusing either on the role of con- 
vective instabilities in the delayed, neutrino-driven explosion 
mechanism within about the first second of evolution (M ez- 
zacappa et al. 1998; Panka & Muller 199^ 

Burrows, Hayes, 

& Fr yxell 1995; [Herant et al. 1994|; M iller/Wilson, & Mayle 
1993), or on the growth of R ayleigh-Taylor instabilitites during 



the late evolutionary stages ([Magataki, Shimizu, & Sato 1998 
Herant & Benz 1992|; [MUller, Fryxell, & Arnett 1991|; Yamada 



& Sato 1991; pachisu et al. 1990[ ). However, multidimensional 
simulations which follow the evolution of the supernova shock 
from its revival due to neutrino heating, until its emergence 
from the stellar surface have not yet been performed. 

In this Letter, we report on preliminary results of high- 
resolution two-dimensional (2D) computations which for the 
first time cover the neutrino-driven initiation of the explosion, 
the accompanying convection and nucleosynthesis as well as 
the Rayleigh-Taylor mixing within the first ~ 300 seconds of 
evolution. 

2. NUMERICAL METHOD AND INITIAL DATA 

We spUt our simulation into two stages. The early evolu- 
tion (f < 1 s) is followed with a version of the HERAKLES 



code (T. Plewa & E. Muller, in preparation) which solves the 
multidimensional hydrodynamic equations using the direct Eu- 
lerian version o f the Piecewise Parabolic Method (C olella & 
Woodward, 1984) augmented by the Consi stent Multifluid Ad- 
vection scheme of Plewa & Muller (1999) in order to guaran- 
tee exact conservation of nuclear species. We have added the 



input physics described in lanka & Muller (1996) (henceforth 



JM96) with the following modifications. General relativistic 
corrections are m ade to the gravitational potential following 
Van Riper (1979). A 14-isotope network is incorporated in or- 
der to compute the explosive nucleosynthesis. It includes the 
13 a-nuclei from ''He to ^''Ni and a representative tracer nu- 
cleus which is used to monitor the distribution of the neutrino- 
heated, neutron-rich material and to replace the ^^Ni produc- 
tion when Kg dro ps below ^ 0.49 (cf. T hielemann, Nomoto, & 
Hashimoto, 1996). Our initial data are taken from the 15 M© 



progenitor model of Woosley, Pinto, & Ensman (1988 ) which 
was collapsed by Bruenn (1993). The model is mapped to a 
2D grid consisting of 400 radial zones (3.17 x 10* cm < r < 
1.7 X lO'^cm), and 180 angular zones (0 < 6* < tt; cf. JM96 
for details). A random initial seed perturbation is added to the 
velocity field with a modulus of 10"^ of the (radial) velocity of 
the post-collapse model. The computations begin 20 ms after 
core bounce and are continued until 885 ms when the explosion 
energy has saturated at 1 .48 x 10^' erg (this value has still to be 
corrected for the binding energy of the outer envelope). We will 
henceforth refer to this calculation as our "explosion model". 

The subsequent shock propagation through the stellar enve- 
lope and the growth of Rayleigh-Taylor instabilities is followed 
with the AMRA Adaptive Mesh Refinement (AMR) code (T. 
Plewa & E. Muller, in preparation). Neutrino physics and grav- 
ity are not included in the AMR calculations. Both do not 
influence the propagation of the shock during late evolution- 
ary stages. However, gravity is important for determining the 
amount of fallback, a problem which is outside the scope of 
the present study. The equation of state takes into account con- 
tributions from photons, non-degenerate electrons, e''"e"-pairs, 
'H, and the nuclei included in the reaction network. The AMR 
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calculations are started with the inner and outer boundaries lo- 
cated at rjn = 10^ cm (i.e. inside the hot bubble containing the 
neutrino-driven wind) and rout = 2 x lO'^cm, respectively. No 
further seed perturbations are added. Our maximum resolution 
is equivalent to that of a uniform grid of 3072 x 768 zones. We 
do not include the entire star but expand the radial extent of 
the grid by a factor of 2 to 4 whenever the supernova shock is 
approaching the outer grid boundary, which is moved from its 
initial value out to rout = 1.1 x 10'^ cm at f = 300 s. Reflecting 
boundary conditions are used at 9 = and d = TT and free outflow 
is allowed across the inner and outer radial boundaries. 

3. RESULTS 

The general features of our explosion model are compara- 
ble to the models of JM96. The most important difference 
is caused by our use of general relativistic corrections to the 
gravitational potential. Since the shock has to overcome a 
deeper potential well than in the Newtonian case, the lumi- 
nosities (which are prescribed at the inner boundary) required 
to obtain a certain final explosion energy are roughly 20% 
higher than those of JM96. In particular, we have adopted the 
following set of parameters: = 2.8125 x 10^^erg/s,L° = 
2.375 X lO^^ei-g/s, Ay, =0.0875,' A£ = 0.0625 (cf. JM96). The 
neutrino spectra and the temporal decay of the luminosity are 
the same as in JM96. 

For the chosen neutrino luminosities the shock starts to move 
out of the iron core almost immediately. Convection between 
shock and gain radius sets in ~ 30 ms after the start of the sim- 
ulation in form of rising blobs of heated, deleptonized material 
(with *C 0.5) separated by narrow downflows with ~ 0.49. 
The shock reaches the Fe/Si interface at r = 1.4 x 10*^ cm after 
~ 100 ms. At this time the shocked material is still in nuclear 
statistical equilibrium and is composed mainly of a-particles 
and nucleons. When the temperature right behind the shock 
drops below ~ 7 x 10'^ K, ^''Ni starts to form in a narrow shell. 
During the ongoing expansion and cooling, nickel is also syn- 
thesized in the convective region. However, this synthesis pro- 
ceeds exclusively in the narrow downflows which separate the 
rising bubbles and have a sufficiently high electron fraction Y^. 
The convective anisotropics therefore lead to a highly inhomo- 
geneous nickel distribution (Fig. |l|). 

Freezeout of complete silicon burning occurs at f « 250 ms, 
and convection ceases at f w 400 ms, at which time the flow 
pattern becomes frozen in. Subsequently, the entire post-shock 
region expands nearly uniformly. The post-shock temperature 
drops below 2.8 x 10*^ K at f = 495 ms, when the shock is about 
to cross the Si/O interface. Thus, our model shows only moder- 
ate oxygen burning (due to a non-vanishing oxygen abundance 
in the silicon shell), and negligible neon and carbon burning. 
Th is is caused by the spe cific structure of the progenitor model 
of Woosley et al. (1988) and may change when different (e.g. 
more massive) progenitors are used. In total, 0.052 Mq of ^^Ni 
are produced, while O.IOM0 of material are synthesized at con- 
ditions with Fe < 0.49 and end up as neutron-rich nuclei. Fig. [l| 
shows the ^^Ni distribution at f = 885 ms. The density ratio be- 
tween the dense regions which contain the nickel and the low- 
density, deleptonized material in the bubbles deeper inside is 
-2.5. 

During the next seconds the shock detaches from the for- 
merly convective shell that carries the products of explosive nu- 
cleosynthesis and crosses the C-nO/He-interface, initially accel- 
erating and subsequently decelerating due to the varying den- 



sity gradient. Twenty seconds after core bounce this unsteady 
propagation speed of the shock has led to a strong compression 
of the metal-containing shell. At the inner boundary of the shell 
pressure waves have steepened into a reverse shock. Rayleigh- 
Taylor instabilities start to grow at the Ni-nSi/O- and C-nO/He- 
interfaces and are fully developed and interact with each other 
at f = 50 s. Nickel and silicon are dragged upward into the he- 
lium shell in rising mushrooms on angular scales from 1° to 
about 5°, whereas helium is mixed inward in bubbles. Oxy- 
gen and carbon, located in intermediate layers of the progeni- 
tor, are swept outward as well as inward in rising and sinking 
flows. At f = 300 s the densities and radial velocities between 
the dense mushrooms and clumps and the ambient medium dif- 
fer by factors up to 5 and 1.3, respectively. As Fig. || shows, the 
fastest mushrooms have already propagated out to more than 
half the radius of the He core. The excessive outflows along 
the symmetry axis in Fig. ^ are a numerical artifact caused by 
the singularity of the polar axis in spherical coordinates. In or- 
der not to introduce large errors, only an angular wedge with 
15° <0< 165° was used in all the analyses. 

The remarkable efficiency of the mixing is clearly visible 
from Fig. || where the original onion-shell structure of the 
progenitor's metal core has disappeared after 300 s. We ob- 
serve large composition gradients between different regions of 
a mushroom and also between different mushrooms. Some of 
them contain ""^Ni mass fractions of more than 70% whereas 
others have only nickel admixtures of 20% or less but show 
high concentrations of silicon and oxygen. Composition con- 
trasts of at least this magnitude have recently been observed in 
different fil aments of the Cas A supernova remnant by H ughes 
et al. (2000). We find that also more than O.O4M0 of matter 
consisting of neutron-rich nuclei are mixed into the ^^Ni and 
^^Si clumps. 

The reverse shock substantially decelerates the bulk of the 
nickel from - 15000 km/s at t = 885ms to 3200 - 4500 km/s 
at f = 50 s, at which time the maximum velocities v™" have 
dropped to 5800 km/s. After another 50 s, however, the clumps 
start to move essentially ballistically through the helium core 
and only a slight drop of v™" from 5000 km/s at f = 100 s to 
4700 km/s at f = 300 s occurs, when most of the ^*Ni has veloc- 
ities below 3000 km/s (Fig.Q). 

We have recently accomplished to follow the subsequent evo- 
lution of our model up to 16 000 s after core bounce. Similar 
to the situation at the C-nO/He interface, the supernova shock 
leaves behind a dense (Rayleigh-Taylor unstable) shell at the 
He/H interface. While the entire shell is rapidly slowed down, 
a second reverse shock forms at its inner boundary (Fig. ^). Our 
high-resolution simulations reveal a potentially severe problem 
for the mixing of heavy elements into the hydrogen envelope 
of Type II supernovae like SN 1987 A. We find that the fast 
nickel containing clumps, after having penetrated through this 
reverse shock, dissipate a large fraction of their kinetic energy 
in bow shocks created by their supersonic motion through the 
shell medium. This leads to their deceleration to — 2000 km/s 
in our calculations, a negative effect on the clump propagation 
which has not been discussed previously. A detailed analysis of 
this result will be subject of a forthcoming publication. 

4. CONCLUSIONS 

Our calculations provide a natural explanation for the mixing 
of ^^Ni into the helium layer, which is required to reproduce the 
light curves of Type lb supernovae ( [Shigeyama et al., 199C| ). In 
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this sense they justify the rather large seed pert urbations (> 5%) 
which were imposed by Hachisu et al. (1994) on the radial ve- 
locity given from spherically symmetric models of exploding 
helium cores at about 10 s after shock formation. In addition, 
they may have interesting implications for the modeling of Type 
lb spectra (Woosley, pri vate communication; comp are Fig. 0of 



this work with Fig. 7 in Woosley & Eastman 1997) 



Our simulations suggest that ballistically moving, metal-rich 
clumps with velocities up to more than ^ 4000 km/s are ejected 
during the explosion of Type lb (and Ic) supernovae. In case of 
Type II supernovae, however, the dense shell left behind by the 
shock passing the boundary between helium core and hydro- 
gen envelope, causes a substantial deceleration of the clumps. 
The high iron velocities observed in SN 1987 A, for example, 
can therefore not be accounted for by our models. This prob- 
lem would be reduced if the density profile at the He/H interface 
were smoother, leading to less strong variations of the shock ve- 
locity. Could this direct to a possible common envelope phase 
or m erger history of the progenitor star (^Hille brandt & Meyer 
1989; podsiadlowski. Joss, & Rappaport 1990| )? Alternatively, 
the high iron velocities in S N 1987 A could requi re additional 
energy input from Ni decay (Herant & Benz, 1992) or could im- 
ply a large global anisotropy of the explosion, e.g. associated 
with j ets emerging from the collapse of a rapidly rotating stellar 
core (Fryer & Heger, 1999). Though final conclusions require 
three-dimensional calculations, our simulations indicate that all 
computations of Rayleigh-Taylor mixing in Type II supernovae 
carried ou t so far (including the case of SN 1987 A; H erant & 
Benz 1992) have been started from overly simplified initial con- 
ditions since they have neglected clump formation within the 
first minutes of the explosion. 

The mixing of nucleosynthesis products as seen in our mod- 
els poses another potential problem. A significant fraction 
(0.04Mq of a total of 0. 1 Mq) of the neutron-rich nuclei which 
are synthesized in regions with < 0.49 is dragged outward 
by the nickel containing mushrooms and clumps. At least this 



amount will most likely be ejected in the explosion. A detailed 
analysis of the nuclear composition is necessary to tell whether 
this is in conflict with limits from Galactic chemical evolution 
models, which allow for at most lO'^'-lO'^M© of neutron-rich 
nucl ei to be ejected per supernova event (e.g., H erant et al. 
1994; [Thielemann et al. 1996| ). Later fallback will not solve this 
potential problem: How could it disentangle the clumpy nickel 
ejecta from their undesirable pollution? A better knowledge of 
the luminosities and spectra of electron neutrinos and antineu- 
trinos emitted from the nascent neutron star is therefore needed 
to perform more reliable calculations of the neutronization of 
the neutrino-heated ejecta. 

During the computations we became aware of oscillations 
with angle in parts of the postshock flow (Figs, [l] and ||). These 
are caused by the "odd-even decoupling" p henomenon associ- 
ated with grid-aligned shocks (Quirk, 1994). As a consequence, 
the maximum nickel velocities, v™", obtained in our AMR cal- 
culations have probably been overestimated by ^25% because 
the growth of some of the mushrooms was influenced by the 
perturbations induced by this numerical defect. The main re- 
sults of our study, however, are not affected. We note that a 
large number of supernova calculations performed with codes 
based either on t he direct Eulerian (cf. Figs. 22 and 24 in Bur- 
rows et al. 1995, Fig. 20 in JM96), or the Lagrangean with 
rema p formulation of the PPM scheme (M ezzacappa et al., 
1998) seem to be affected by this numerical flaw. We defer a 
detailed analysis of this problem to a forthcoming publication. 

We thank S. Woosley for profiles of the progenitor star and 
S. Bruenn for his post-bounce core model. We acknowledge 
support by P. Cieciel^g and R. Walder concerning visualization 
and by the crew of the Rechenzentrum Garching where the sim- 
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computers. TP was supported by grant 2.P03D.004.13 from the 
Polish Committee for Scientific Research, HTJ by DFG grant 
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Fig. 1 . — Distribution of the '*Ni mass fraction at f = 885 ms. At this time the shock (not visible) has ahnost reached the outer boundary of the computational 
domain (blue) at 1.7 x lO' cm. 
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Fig. 2. — Density distribution (logjQp [gem ^] ) 300 s after eore bounce in the inner ^ 3 X lO" cm of the star. The supernova shock (outermost blue discontinuity) 
is located inside the hydrogen envelope at r = 2.7 X lO" cm. A dense shell (visible as a red ring) has formed behind the shock. Its outer boundary coincides with 
the He/H interface while its inner boundary is in the process of steepening into a reverse shock. 
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Fig. 3. — Initial composition exterior to tlie iron core (top) and composition 300 s after core bounce (bottom). The CVO core of the progenitor lias been 
completely shredded by the Rayleigh-Taylor instability. '^C, '*0, ^*Si and the newly synthesized '*Ni have been mixed beyond the inner half of the helium core 
(Mne-core = 4.2Mq) to a mass coordinate of 2.45 Mq. Inward mixing of '^He has led to an increased helium concentration interior to 1.75 Mq. Inside 1.4 Mq the 
a-rich freezeout in the neutrino-heated ejecta yields a contribution of '^He. A cone with opening angle of 15 degrees around the polar axis was excluded from the 
analysis (see text). 
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Fig. 4. — Distribution of the '''Ni mass versus velocity att = 50 s and t = 300 s, respectively. A cone with opening angle of 15 degrees around the polar axis was 
excluded from the analysis (see text). 
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